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ABSTRACT 

Recently, the gravitational instability (GI) model of giant planet and brown dwarf formation 
has been revisited and recast into what is often referred to as the "tidal downsizing" hy- 
pothesis. The fragmentation of self-gravitating protostellar discs into gravitationally bound 
embryos - with masses of a few to tens of Jupiter masses, at semi major axes above 30 - 40 
au - is followed by a combination of grain sedimentation inside the embryo, radial migration 
towards the central star and tidal disruption of the embryo' s upper layers. The properties of the 
resultant object depends sensitively on the timescales upon which each process occurs. There- 
fore, GI followed by tidal downsizing can theoretically produce objects spanning a large mass 
range, from terrestrial planets to giant planets and brown dwarfs. Whether such objects can 
be formed in practice, and what proportions of the observed population they would represent, 
requires a more involved statistical analysis. 

We present a simple population synthesis model of star and planet formation via GI and 
tidal downsizing. We couple a semi-analytic model of protostellar disc evolution to analytic 
calculations of fragmentation, initial embryo mass, grain growth and sedimentation, embryo 
migration and tidal disruption. While there are key pieces of physics yet to be incorporated, it 
represents a first step towards a mature statistical model of GI and tidal downsizing as a mode 
of star and planet formation. 

We show results from four runs of the population synthesis model, varying the opac- 
ity law and the strength of migration, as well as investigating the effect of disc truncation 
during the fragmentation process. We find that a large fraction of disc fragments are com- 
pletely destroyed by tidal disruption (typically forty percent of the initial population). The 
tidal downsizing process tends to prohibit low mass embryos reaching small semimajor axis. 
The majority of surviving objects are brown dwarfs without solid cores of any kind. Around 
forty percent of surviving objects form solid cores of order 5-10 Earth masses, and of this 
group a few do migrate to distances amenable to current exoplanet observations. Over a mil- 
lion disc fragments were simulated in this work, and only one resulted in the formation of a 
terrestrial planet (i.e. with a core mass of a few Earth masses and no gaseous envelope). These 
early results suggest that GI followed by tidal downsizing is not the principal mode of planet 
formation, but remains an excellent means of forming gas giant planets, brown dwarfs and 
low mass stars at large semimajor axes. 

Key words: 

planets and satelUtes:formation, stars; formation, accretioniaccretion discs, methods: numeri- 
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1 INTRODUCTION 

Traditionally, there have been two distinct models or paradigms of 
planet for mation in proto s tellar discs. In t he co re accretion (CA) 
paradigm JPollackl 119961 : iHubickvi et alj 120050 . planetary mass 
bodies are assembled in a bottom-up process. The protostellar disc 
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initially contains a population of interstellar grains of size ~ l/iin, 
which coagulate into larger particles. As the grains grow past a 
few cm in size, they settle to the midplane, aggregating into plan- 
etesimals a few km in size. The planetesimals can grow via colli- 
sions into solid protoplanetary cores, which can then accrete a thin 
gaseous envelope. 

The fate of this solid core then depends on the mass it can 
reach - if the core exceeds a critical core mass of around 10 
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Earth masses ( lMizundll98 0). the envelope accretion process is no 
longer thermally regulated, resulting in runaway growth of the at- 
mosphere, producing gas and ice giants. Below this critical core 
mass, the gas accretion remains slow, and teiTestrial planets can be 
formed. 

Core accretion is certainly the most accepted model of planet 
formation, and it is the most successful at explaining the major- 
ity of the observed exoplanet population. Despite this, predictions 
of deserts in the exoplanet distribution at short periods and low 
mass have not been borne out by observatio ns, and the existence 
of highly compact systems such as Kepler- 11 JLissauer et alj201 ih 
are not easily explained. 

Questions also remain regarding certain stages of the core 
growth process. There are barriers to grain growth at various size 
scales. The most famous is the "metre baiTier", which is in reality a 
mixture of inhibiting factors that take effect for grain sizes as small 
as millimetres depending on the local disc properties. As grains 
increase in size, their relative velo cities increase, and th e probabil- 
ity of continued growth decreases ( IBlum & Wurml2008h . Also, gas 
drag works most efficiently on metre sized bodies, caus ing rapid 
radial drift and loss of solid m aterial onto the central star IWhippld 
I1973I : IWeidenschillindl 1977h . More recently, light has been shed 
on the "bouncing barrier", which occurs at slightly lower rela- 
tive velocities, and als o inhibits grain growth ( IGuttler et al.ll20ld : 
IWindmarket"ai]|2012l) . 

Conversely , in the gra vitational instability (GI) paradigm 
JCamerod Il978l : IBossI 1 19970 . planetary-mass objects are assem- 
bled in a top-down fashion through fragmentation of the gaseous 
protostellar disc. This occurs during a relatively brief phase of the 
disc's life where it is sufficiently massive to be self-gravitating, and 
can become gravitationally unstab le. Gravitation ally unstable discs 
must satisfy the Toomre criterion dToomrdl 19641) : 



Q 



TvGT, 



< 1.5-1.7, 



(1) 



where Q is the Toomre parameter, Cs is the local sound speed of 
the gas, Ke is the local epicyclic frequency (equal to the angular 
velocity fl in Keplerian discs) and E is the disc surface density. 
This is a linear stability criterion, and for axisymmetric perturba- 
tions Q < 1 for instability. Equation ([T) refers to non-axisymmetric 
perturbations, where the critical Q valu e has been established by 
simulations (see e.g. lDurisen et al.l2007l for a review). 

The onset of gravitational instability produces trailing spiral 
density waves which heat the disc, increasing Q. The density per- 
turbations produced by the spiral waves eventually decay, allow- 
ing the disc to cool and become unstable once more. As new spi- 
ral structure is generated by the instability, this cycle of growth 
and decay introduc es turbulence in to the gas, usually referred to as 
gravito-turbulence ( Gammidl200lh . If the disc's angular momen- 
tum transport is locally determined, this turbulence can be rep- 
resented by a "pseudo-viscosity" v ([Balbus & Papaloizod 1 19991 : 
iLodato & Ricdl2004 : Forgan et alJuOl 11) , and as such is amenable 
to the lShakura & Sunvaed ( Il973h a-parametrisation: 



u = aCsH, 



(2) 



where H is the disc scale height. This description greatly simplifies 
the disc evolution to a viscous form (although typically non-linear 
and only soluble numerically) such that, if we assume that the disc 
is axisymmetric, then the surface density E(r, t) evolves according 
to JLvnden-Bell & Pringldl 1 9741 : iPringlell 1 98 lb 
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where Ewind represents mass loss from the disc due to a disc wind. 
The disc can respond to the heating caused by the instability via 
radiative cooling, reducing the local sound speed and decreasing 
Q. As self-gravitating discs will only exist in the early stages of 
star formation, we should also expect a mass loading term in these 
equations which can have significant consequences for the disc's 
evolution (see e.g. Clarke 2009) . 

As a result of the competition between heating and cooling, 
gravitationally unstable di scs can maintain what is referred to as a 
marginally unstable state jPacz vnskill 1 97 8h . where the shock heat- 
ing and radiative cooling become balanced, Q maintains a self- 
regulated value around the critical value, and the disc becomes 
quasi-steady. If the disc is marginally unstable a nd in local the r- 
modynamic equilibrium, then we can solve for a ( lGanmiidl200lh : 



(4) 



97(7 - l)/3c 



where we have defined /3c as a dimensionless expression of the 
local cooling time tcooi : 



Pc — tcool^'. 



(5) 



For a self-gravitating disc to fragment, this self-regulating mech- 
anism must be broken. The density perturbat ions induced by the 
gravitational ins tability have the behaviour (ICossins et al.l l2009l : 
lRiceetal.ll201lh : 



AEr- 



>cx 



(6) 
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Therefore, to produce a density perturbation of sufficiently high 
amplitude that it can collapse into a bound object, the local cooling 
time must be sufficiently short, or equivalently the local gravita- 
tional stresses must be sufficiently high. 

It is reasonably clear that self-gravitating stresses saturate at 
some value of a, setting a lower limit for the cooling time, below 
which a marginally unstable state can no longer be maintained. Be- 
low this critical cooling time /3crit, the disc is no longer able to 
achieve thermodynamic balance, and fragmentation is expected to 
occur. The critical cooling time is a function of the ratio of spe- 
cific heats (7), but the saturating value of a is fixed at around 0.06 
l lGammidl^0ll : lRice et ai]|2005l) . 

The precise satur ating value of a is currently unclear due 
to convergence issues JMeru & Batd l201 1). Indeed, it may be the 
case that fragmentation proceeds in a stochastic fashion, where the 
power spectrum of density fluctuations is such that a weak pertur- 
bation results in fragmentation, and the traditional fragmentation 
criteria do not delineate fragmenting/non-fragmenting configura- 
tions, but indicate where the probability of fragmentation is close 
to unity. If this is the case, it may be that fragmenta tion can occur 
with low probability at values of Pc a s high as 20-50 jPaardekoopen 
l2012l : lHopkins & Christiansenl2013h ). This debate is yet to be fully 
resolved, but there are indicat ions that carefu l consideration of the 
effects of artificial viscosity JLodato &Clarkdl2pil[) and a more 
appropriate implementation of radiative cooling dRice et al.ll2012h 
can help to address this problem, and that the saturating value of a 
remains of order 0.1. As yet, true stochastic fragmentation has not 
been observed in global disc simulations. Due to our current igno- 
rance of how its properties may be incorporated into this work, we 
will not include it in this analysis. 

Regardless, these arguments have to some extent excl uded 
the effects of mas s infall, thermal hist ory (Clar ke et al.ll2007 ) and 
stellar irradiation dKratter & Murrav-Clav 2011) . Later work has 
generalised this fragmentation criterion in terms of critical scales 
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such as the Jeans Mass ('Forga n & Ric3l201ll . l2013h or the Hill Ra- 
dius l IRogers & Wadsley 2012) , which agree well wit h results de- 
rived from direct investigation from the spiral structure (IBolev et alj 

l2oToh . 

Requiring the Q-condition for gravitational instability and a 
second criterion for fragmentation restricts the disc radii at which 
fragments may form. Typically, the inner disc does not cool effi- 
ciently enough to maintain low enough Q values, and is therefore 
stable against self-gravity. Also, gravitational stre sses tend to de- 
crease with increasing proximity to the central star ( Armitage et al.l 



I2OOII ; iRice & Armitag3l2009l : Iciarkell2009l : iRice et al.ll201 

Fragmentation is therefore consigned to the oute r disc, 
r adii a bo ve 30 or 40 au (IRafikovl 
'2005": 



and typically to 



Matzner & Levin 
Bolevetal. 2006 



Whitworth & Stamatellos 



2005 



2006; 



Stamatellos & Whitworth 2008; Fo rgaii et all 
2009; Clarke 2009; jKratter et al.l I2OI0I; IVorobvov & Basul I2OIOI; 






Forgan & Rica .20111) . Estimat es for the fragment mass tend to 



place a lower limit of 3-5 Mjup JKratter et alj2010l ; lForgan & Ricd 
uOllh . These two results alone would appear to be enough to dis- 
count the gravitational instability theory as a means of forming 
lower mass planets at disc radii inside 30 au. Equally, core 
accretion has difficulty explaining the existence of massive planets 
with semi-major axes above 30 au, and so it may be the case that 
planet formation proceeds via both core accretion and gravitational 
instability (Boley 2009). 

However, the recent reformulation of the gravitational insta- 
bility theorjMnto what isconunonl^;_refa downsiz- 
ing" JBolev et al.ll2oTol. I2OI ll ; lNavakshinll20Tli. l2010al lbb has re- 
vived the possibility that the GI formation mode could produce low 
mass planets at low semimajor axis. The key to this revival is the 
study of the subsequent evolution of the disc fragments into plan- 
etary embryos. Three key physical processes shape the fate of the 
embryo: 

(i) Dust growth and sedimentation. As grains grow inside the 
embryo, the sedimentation velocity of the grains increases, draw- 
ing them towards the pressure maximum at the embryo's centre. If 
this process is efficient, the gravitational field at the centre of the 
embryo will become dominated by solids. The solids can become 
self-gravitating, potentially even gravitationally bound, forming a 
solid core. As the temperature at the centre of the embryo increases, 
the grains can eventually be vapourised, ending the growth of the 
core. 

(ii) Radial migration. Differential torques induced by the disc 
on the embryo can result in inward radial migration, in much the 
same way as protoplanets in the core accretion paradigm. 

(iii) Tidal disruption. As an embryo migrates inwards, its phys- 
ical radius will be comparable to the Hill Radius. As this occurs, 
the upper layers will become unbound, and hence stripped from the 
embryo. 

By comparing the time scales on which these three processes oc- 
cur, iNavakshinI yOlOaj) gives several possible outcomes for the 
embryo. In some cases, it is conceivable that a solid core of 
less than ten Earth masses can be formed inside the embryo, 
which has its upper layers completely disrupted, leaving a ter- 
restrial protoplanet. Equally, inefficient core formation and strong 
tidal disruption could result in the entire embryo being destroyed 
and accreted by the sta r , giving rise to outburst behaviour (e.g. 
Vorobyov & B asul |2005|; iBolevet al.l I2OI0I ; IPunham & Vorobvovl 
20I2l ; INavakshin & Lodatdl201 2). Intermediate cases produce gi- 
ant planets both with and without cores depending on the sedimen- 
tation timescale. If the core is not fully gravitationally bound be- 



fore the embryo is disrupted, then this could potentially produce 
belts of planetesimals (iNayakshin & Cha 2012 ), or generally de- 
liver processed solid material that can affect the overall disc chem- 
istry l Bolevetal.l2oToll201lh . 

This production of planets of varying types and masses, which 
circumvents the need for planetesimal growth without precluding 
the existence of asteroid belts, makes the tidal downsizing hypoth- 
esis an attractive one for planet formation theorists. However, it 
remains unclear what distribution of object mass and semi-major 
axis should be expected from self-gravitating disc fragmentation. 
A common criticism of gravitational instability theories as a whole 
is the lack of population synthesis models, which have been im- 
plemented forJhe_£ore_accretimjnodelvwt^^ great deal of suc- 
cess (e.g. [I da & Lin' l 2004 lAlibert et all l2005l ; 1lda & Liril I2OO8I; 
iMordasini e t al. 2012). 

In this paper, we wish to remedy this by presenting a simple 
population synthesis model for disc fragmentation followed by tidal 
downsizing. This model is still missing some physical processes, 
and simplifies others, but it gives the broad strokes of what stellar 
and planetary systems we can expect to see produced by this the- 
ory. As such, we present it as the first stage in a continuing series of 
work (much as was done with the population synthesis models of 
core a ccretion). We couple the equations of Navak shinI (l2010al lbl 
uOIlal) with semi-analytic m odels of disc evolution incorporat- 
ing X -Ray photoevaporation ( IRice & Annitagell2009l ; lOwen et al.l 
l20Ilh . and evolve a large ensemble of planetary embryos formed in 
a variety of self- gravitating protostellar discs. While uncertainties 
remain about the underlying distributions of disc mass and radius 
at early times, a uniform sampling of these parameters allows us 
to explore the exoplanet parameter space in which we can expect 
the GI formation mode to leave its signature. This will allow the 
statistical properties of the two principal planet formation theories, 
core accretion and GI, to be compared much more directly than has 
been possible until now. 

The paper is structured as follows: in section |2]we discuss the 
implementation of the equations that describe the sedimentation, 
migration and disruption of embryos, as well as the disc models 
used in this work. In section |3] we display the resulting statistical 
distributions as a function of a few remaining free parameters. In 
section |4] we discuss the implications of these models for the tidal 
downsizing hypothesis, critically assess the current weaknesses of 
our approach and advocate avenues for further work. Finally, in 
section|5]we summarise the work. 



2 METHOD 

2.1 A simple semi-analytic model of tidal downsizing 

In this section we summarise the equations governing embryo for- 
mation, core formation and the subsequent migration/disruption 
processes that can occur. Throughout, we rely heavily on the equa- 
tions given by Nayaks hin in his three pa pe rs outli ning the tidal 
downsizing hypothesis ( lNavakshinll2010a| b |. 201 laf). We also us e 
semi-analytic disc models as described in lRice & Armitagel ( l2009f) . 
and determine fragmentation a nd fragment masses us ing the Jeans 
mass formalism as described in lForgan & Ricd (1201 11) . 
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2. 1. 1 Conditions for Disc Fragmentation 

To firstly assess where fragmentation can occur, we calculate the 
Toomre Q criterion at all disc radii, as well as the Tj criterion 
outlined in lForgan & Ricd ( 120 111) , where 



Mj 
Mj 



n. 



(7) 



and Mj is the local Jeans mass inside a spiral perturbation. This 
condition allows us to determine at what rate the Jeans mass is 
changing. Regimes where the Jeans mass is decreasing rapidly will 
be prone to fragmentation. If we assume the disc is marginally un- 
stable and in a steady state, then 



P _ ; 0/9/ 1 , 9Q7(7-1) 
r,-(3/2(~- + 



(8) 



If the disc is in thermal equilibrium, then Tj — )■ 00, and 
the disc will not fragment. If /3c is sufficiently low that thermal 
equilibrium would demand that a would exceed the value asat 
at which gravito-turbulent transport is expected to saturate, i.e. 
Of > Qsat ~ 0.1, we set a — cvsat and Fj becomes negative. As 
Pc decreases, Tj remains negative, and its magnitude decreases. 
Once Fj is sufficiently small and negative, we can see that the lo- 
cal Jeans mass is decreasing sufficiently quickly that a local spiral 
density perturbation will soon exceed the Jeans mass, and begin the 
fragmentation process. 

We therefore demand that (in disc regions where Q ^ 1.5) the 
local Jeans mass must be decreasing at a rate such that 



5 < r./ < 



(9) 



for fragmentation to occur. The exact critical values for Tj are cur- 
rently unclear, but the above pr oduce results consistent with con- 
ventional cooling time criteria (IForgan & Ricell201 iL 120131) . The 
value of Fj sets a timescale for fragmentation rather than spec- 
ifying whether or not it will occur, if Fj < —5, then in theory 
fragmentation remains possible, but unlikely, as disc evolution will 
typically act to prevent this. 

The advantage of using the Jeans mass criterion is that 
it immediately iden tifies an initial fragment mass M. From 
iForgan & Ricj ( 1201 ih . this is: 



M = Mj 



3G (1 + Ap)' 



(10) 



While disc fragmentation is an inherently three dimensional pro- 
cess, the initial conditions provided by the disc are roughly two- 
dimensional. If the surface density perturbation induced by a spi- 
ral density wave AE/E is of order unity, the perturbations will 
become non-linear, and fragmentation can be expected to occur 
(this statement is true despite the current debate regarding stochas- 
tic fragmentation). In other words, at the point of fragmentation, 
a 2D measure is sufficient to identify when fragmentation occurs, 
but the details of the fragment's subsequent evolution are 3D. As 
the discs we use in the modelling done here are axisymmetric, sur- 
face density is in fact a ID measure, and for these purposes is suf- 
ficient to control/govern the fragmentation process. As with much 
semi-analytic modelling, using ID measures as is done here is a 
necessary and appropriate sacrifice for the sake of computational 

expediency. 

Using the empirical relation for AE/E derived bv lRice et al.l 

( 1201 ih , we can give the following form for M: 



M 



3G (1 + 4.47^/S) ' 



(11) 



Fragments are formed wherever both the Q and Jeans mass criteria 
specified above are satisfied. We might expect fragments to foim at 
spacings of appr oximately 2 Hil l radii, so that feeding zones do not 
strongly overlap (ILissaueiil987h . To reflect our current uncertainty 
in this matter, the precise value of the spacing is uniformly dis- 
tributed between 1.5 and 3 Hill radii. This also allows us to probe a 
larger variety of initial fragment properties from what is a limited 
set of initial disc models. 

This is sufficient to create the initial embryo population. We 
must now investigate the evolution of the embryos individually. 

2.1.2 Initial properties of the fragment/embryo 

By appealing to the analogy of t he opacity limit for fragmentation 
in star formation (Rees 1976 ; Low & Lvnden-Bell 1976). it can be 
seen that fragments are initi ally akin to the "first co res" (FCs) of the 
star formation process (e.g. lMasunaga et al.ll998n . We will refer to 
these as FCs to avoid confusion with the later solid cores that may 
develop inside these embryos. 

If these fragments are indeed FCs, it is appropriate to model 
them as poly tropic spheres. The properties of the poly trope depend 
on the initial temperature (taken usually to be ~ IQK) and the crit- 
ical density pad at which the gas switches between isothermal and 
adiabatic behavi our. The se quantities const rain the initial ra dius of 
the FC cMasunaga et alj|T998 - Masunaga & Inutsukall 19991) . If we 
follow the route of lNavakshinI ( l2010al) and define the opacity law 

and the following scaled quantities: 

M 
mi 



Ti 



O.OIM0 
T 



WK 

aoT' 



we find the FC radius R as (see lNavakshinll2010al for more): 



R= 17.5a.um~''-'T^ 



l/3^(l+4p«)/9 4/9 



(13) 



(14) 



The temperature of the FC will be close to the virial temperature: 



1 vir — 



GMfi 



= UQKm'l'-'T; 



4/3^-{l+4p„)/9 -4/9 



(15) 



3kBR 

The luminosity of the FC is a non-trivial function of the opacity 
power law index p^, and becomes somewhat cumbersome to write 
down in general terms. Hence, Nayakshin evaluates expressions for 
Pk — 1,2, and as such gives the equivalent expressions for the 
cooling time for both opacity laws, assuming the binding energy of 
a sphere: 



J- oon ^2/3rrn-4/3 1/9 

tcooi.o = 5700yrm?r-'^/\;'/^ 



(16) 
(17) 



The "0" in the subscript indicates that these are initial values for the 
cooling time. Nayakshin shows that as the FC contracts, the cooling 
time increases: 



icool(i) — icool.O + (1 + Pk.) t. 

As a result, the temperature evolves according to: 

t 



T{t) = r(o) 



i + (i+P«) 



(18) 



(19) 
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The evolution of the temperature will be crucial to the survival of 
grains inside the embryo, as will be discussed in the next section. 
Finally, we can express the contraction of the embryo as 



i?(i) = 



R{0) 



[1 + 2^^1 

L tcool.oj 



1/2- 



(20) 



These analytic expressions s hould be cont r asted with the more de- 
tailed calculations made by iHelled et alj <200Q1 a nd subsequent 
works jHelled & Schubertll2008l : lHelled et al.ll2008l) which numer- 
ically integrate the equations of stellar evolution, with a spherically 
symme tric embryo as initial conditions. Figure 1 of Nayak shinl 
( 1201 Obi) shows typical temperature evolution tracks using the above 
equations for FCs between 3 and 30 Mjup, and these com- 
pare well w i th the temperature evolution shown in Figure 2 of 
iHelled et alj (200a). However, we should note that the contraction 
timescales shown here are much longer than those given by calcu- 
lations where the opacity is al lowed to vary due to grain evolution 
( IHelled & Bodenheimenl201 ll) . and this has important implications 
for future work (see Discussion). 



2.1.3 Grain growth and sedimentation 

If there has not been significant grain growth or enrichment before 
disc fragmentation, the protostellar disc grains (of no more than a 
few microns in size) will be extremely slow to sediment to the cen- 
tre of the embryo. Therefore, an initial phase of grain growth must 
be completed before the embryo can form a core on a reasonable 
timescale. This grain growth ti mescale f gr is initially defined b y the 
Brownian motion of the grains dDuUemond & Dominikll2005l) . but 
quickly becomes dominated by differential settling, where (at rel- 
atively small grain sizes) the larger grains experience greater sedi- 
mentation velocity. 

Initially, the grains are small enough for their motion to be 
tightly coupled to that of the gas. The grain growth process is phys- 
ically very similar to that of the core accretion paradigm, with 
an added gravitational acceleration towards the embryo's centre, 
which for a constant density embryo is: 



(^3r{r) = -4/37rGpr, 



(21) 



where r = [0, i?] is the grain's radial position. The drag force due 
to the embryo's gas can be expressed as a combination of Stokes 
and Ep stein regimes ( IWeidenschilling|l977l : lBossl 19981 : iNavakshinI 
l2010al) : 



dVa 

~dt 



pCs 



X 



T,s X + s 



[vs — v) — A/SnGpr, 



(22) 



where A is the mean free path of the grain, s is the grain size, and 
Vg and V are the velocities of the grains and gas respectively. We 
have defined Es = pbuikS, where the bulk density of the grains 
Pbuik ~ 1 gcm~'^. Grains will quickly reach a terminal velocity 
iiscd, which can be found by solving for dv^/dt = 0, 



^scd — 



~47rGEs7?A + s 



3cs 



A 



(23) 



Differential sedimentation allows the larger grains to sweep up the 
smaller grains, and results in a mass growth rate: 



dm 2r 

— =nsf,pv^. 



(24) 



where s is the grain size, fg is the mass fraction in grains, and p 
is the local gas density (i.e. in this initially well-mixed phase, the 



volume density of grains pg = fgp). Converting this into a growth 
rate for s gives: 



ds 
di 



faP. 



-r = -; — "sod, 



and so we can define a grain growth e-folding timescale: 



_ _s_ _ 3ca 

■= - ^ - .fgpGR- 



(25) 



(26) 



Given an initial grain size (so) and final grain size (Scrit) before 
sedimentation can become significant, then the growth timescale is 



— 3^3 
""" ^ irfgpGR 



In 



So 



(27) 



Once the grains have grown to Scrit, the entire body of grains 
(which we refer to as the "grain sphere") will begin to sediment, 
quickly accelerating to critical velocity. As the sedimentation ve- 
locity is a function of R, the sedimentation time is independent of 
radius, i.e. all shells of sedimenting dust will reach the centre of the 
embryo at the same time (if the embryo is uniform). The sedimen- 
tation time is 

R 3cs X 



f-scd 



WsGd 



47rGE^ X + s 



(28) 



The s edimentation of grain s can be limited by destructive colli- 
sions ( IBlum & WurmI 120081) . At velocities of several metres per 
second, collisions will result in fragmentation rather than growth. 
This limits the maximum value of tised, an d as such the sedimen- 
tation timescale has a lower limit (INavak shin 2010b). Generally 
speaking, these equations show similar dependencies on tempera- 
ture and dust t o gas ratio as those of Helled' s calculations (see e.g. 
Appendix 1 of IHelled & Bodenheimenl201ll) , but we note that the 
physics we employ here is much simpler. 

In particular, note that our equations describe the bulk mo- 
tion of the grains, and do not incorporate the changing radial pro- 
file or temperature of the embryo. The evolution of the clump's 
photosphere will have consequences that these equations cannot 
yet model. Radiative heating of the embryo can have significant 
consequences for grain settlin g and destruction (iHelled et al.j2008l : 
IHelled & Bodenheimeiir201 ih . 



2.1.4 Turbulence as a suppressor of settling 

How should we select Scrit in equation (|27|l? We have neglected 
turbulen ce in this analysis, which wil l act to prevent gravitational 
settling JFromang & Papaloizoul200a) . We can approximate the ef- 
fects of turbulent mixing using a ID diffusion equation in spher- 
ical geometry, with the d iffusion coefficient selected using the 
IShakura & SunvaevI ( 119731) viscosity prescription: 



D = aCsR. 



(29) 



If the gas is initially a constant density sphere, the dust sphere set- 
tles to an equilibrium solution 



(30) 



/5s = Ps 

where 
Hd = \ 


oexp 


r R^] 

[ m\ 


•> 


'2DU, 


d- 


If we re-express this as an aspect ratio: 


Hd 


'Sa E A 


1/2 


R 


2 E. 


X + S 





(31) 



(32) 
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For gravitational settling to dominate over the turbulent mixing the 
above aspect ratio should be much less than 1 . For a given turbu- 
lence parameter a, we can assume that Scrit is found by setting 
H d/R = 1 and sol ving for s. Unlike the detai l ed ca lculations 
of iHelled et aTl J2008l) and lHelled & Bodenheime3 J201lh . convec- 
tion is not considered here. They find that embryos are likely to be 
highly convective initially, and this has important consequences for 
slowing the settling process, and returning settled grains back to the 
upper envelope. It is therefore crucial that future models incorpo- 
rate convection (see Discussion for more). 



2.1.5 The ejfect of the contracting embryo 

The grain growth timescale in equation ( 127) tacitly assumes a static 
embryo. However, the grain growth timescale is sufficiently long 
compared to the cooling time that we should expect the embryo to 
cool and contract significantly during this time. This compression 
should accelerate grain growth. The growth timescale scales as 



Cs fT{t)\ 
" pR \t{0)J 



-3/2 



Altering the size growth equation as follows: 



ds 
dt 



l + (l+p„ 



t 



icool.O 



1.5/(l+p^) 



(33) 



(34) 



Integrating this gives the more accurate estimate of the growth 
timescale: 



1+1 



1 I 2.5+p. 



(35) 



where ^ = (1 + Pk)/(2.5 + pn). This revised growth timescale 
can be over an order of magnitude shorter than the static estimate, 
tgr,o- Indeed, numerical experiments indicate that the true timescale 
can be another factor of two shorter still (Navakshin 2010a), and 
as mentioned already, these timescales can be orders of magni- 
tude shorter when opacity evol ution due to grain growth is included 
JHelled & Bodenheimedbouh . 

2.1.6 Core Formation 

The total mass in solids is 

Ms = f,M, (36) 

where fg — 0.02 for solar metallicity. For an embryo of 1 i\fjup, 
this gives a fiducial maximum core mass of around 60Af g . Approx- 
imately two thirds of the solids mass is ices (Navakshin 2010b), so 
melting can reduce this maximum mass by the same factor. 

As the sedimentation timescale is independent of radius, the 
grain sphere, initially of radius Rg (0) (which is equal to the radius 
of the FC), will begin contracting according to 



Rg{t) = Rexp 



t 

f'scd 



(37) 



This allows us to identify the timescale fseif on which the density 
of grains exceeds the local gas density, and the grains become self- 
gravitating. This timescale is 



f'sclf — L^scd 



In/; 



(38) 



The factor of 3 emerges from assuming a spherically symmetric ge- 
ometry. Once the dust sphere is self-gravitating, the sedimentation 
velocity is modified so that the gravitational acceleration comes 
only from mass inside the dust sphere: 



^scd 



GMs{Rs) E, A + s 



i?2 



pCs A 



(39) 



The mass interior to the dust sphere remains constant, and therefore 
the differential equation for Rs takes the form 



dRs 
dt 



A_ 
RV 



where A is some constant, and has the solution: 



R{t) = i?(isolf) 



.t-ts 



f'scd 



1/3 



(40) 



(41) 



The contraction remains homologous, but is more rapid than the 
non-self-gravitating stage. The free fall time of the dust sphere, 
i.e. the time required for all grains to collapse to infinite den- 
sities, becomes ap proximately twice the sedimentation timescale 
dNavakshinlbOlOd) . 

In this self-gravitating phase, the grain sphere can still be 
dissipated. If, for example, the outer envelope of gas is removed, 
the resulting pressure gradient in the grain-gas mix can destroy 
the grain sphere. This is in part because the grain sphere is self- 
gravitating, but not gravitationally self-bound. If the grain sphere 
can contract sufficiently to become self-bound, then even the dissi- 
pation of the outer envelope will not unbind the sphere. Navakshinl 
( 1201 lah refers to this stage as the "grain cluster" stage, as an anal- 
ogy to star clusters which survive interstellar gas dispersal. 

If the grain sphere's self-gravity dominates over the local gas 
pressure gradient, then the grain sphere can enter the grain cluster 
phase. This is in essence a Jeans length argument, and the usual 
calculations apply. The pressure gradient is of order 



dP 



pksT 



dR prriHR' 



(42) 



and should be equal to the gravitational force per unit volume for 
hydrostatic equilibrium: 



, GM,nc{R) , , , 

/grav ~ -^ (P + Ps) ■ 



(43) 



Equating these (and assuming p « pg) gives the Jeans length as 



\.j = 



3 fcsT p 



AttG pniH pi 



1/2 



(44) 



Once the grain sphere radius contracts below this length, then the 
final process of core formation can begin, which is presumed to be 
rapid (Boss 1998; Navakshin 201 la), although the process requires 
further study to ascertain this (see Discussion). 



2.1.7 Ice Melting, Grain Vapourisation and Hydrogen 
Dissociation 

A fundamental barrier to solid core formation is the increasing tem- 
perature of the embryo. As the embryo contracts, the central tem- 
perature increases, and the chemistry and phase of the grains will 
change. At temperatures above ~ 200A", the ice mantles of the 
grains will melt, preventing volatiles from becoming part of the 
core. At temperatures above 1400 K, the grains can be vapourised, 
preventing a solid core from forming at all. Embryos of relatively 
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large initial mass will have higher initial temperatures, and as a re- 
sult will be more likely to vapourise their grains before a core can 
form. 

At temperatures above 2000 K, the molecular hydrogen inside 
the embryo begins to dissociate. The energy that must be stored 
in internal degrees of freedom to achieve this depletes the em- 
bryo's pressure support, allowing it to begin a second collapse. In 
star formation, this pre cipitates the formation of the "second core" 
JMasunaga et alj|l998 !). the bound entity which is the protostar it- 
self. Once this has occurred, the first core equations no longer ap- 
ply, and the embryo has attained a new quasi-static configuration 
with a greatly reduced physical radius. 

Wh en objects reach this s tage, we switch to the mass-r adius 
model of iBurrows et all ( 1 19971) (see e.g. i Burrows et allbOllh . We 
deem any object which reaches this stage with a total mass greater 
than ISAfjup to be a "brown dwarf". Objects which reach this 
stage and have a mass less than ISMjup are split into several cat- 
egories depending on their final configuration. Objects which have 
managed to form a core before melting are recorded as "icy core" 
objects; objects which form a core after melting are recorded as 
"rocky core" objects; and objects which do not form a solid core 
before the grains vapourise are simply labelled "no core" objects. 

From equation dlSt . T oc Af^''', which gives a lower mass 
limit for objects which are able to undergo H2 dissociation before 
being tidally disrupted. 



2.1.8 Disc Migration 

While the embryo is evolving internally, the disc is evolving ex- 
ternally, and differential torques will result in planetary migra- 
tion. The current wealth of lit erature on planetary m igration (see 
iBaruteau & Massell |2013| and iKlev & NeJsor] I2OI2I for reviews) 
is predicated on planets represented either as points or as uni- 
form density spheres. The "first core" calculations shown above 
negate both assumptions, as the embiyo radius will typically ex- 
ceed the Hill radius, and as such more sophisticated hydrodynamic 
simulations are required to show the migration of disc fragments 
in a disc which may still be self-gravitating (e.g. ISaruteau et al.l 
I2OIII ; lMichaeletal1l201ll ; Icha & Navakshinl I2OIII ; IBossI|2013|) . 
These simulations have shown that the migration of clumps in 
self-gravitating discs is non-trivial, and typically much more rapid 
than calculations based on non self-gravitating discs would pre- 
dict. Also, the migration direction becomes a function of the disc's 
viscosity, as has been found in other migration simulations (cf. 
lBitschetall2013h . 

This being said, it is useful to us e the standard mi gration 
timescale calculations (as was done by iNavakshirj l2010al) to al- 
low comparison with their work, as well as to keep the complexity 
of the semi-analytic model under control. Indeed, it is as yet un- 
clear if a full semi-analytic description of fragment migration has 
emerged from numerical simulations, and will most likely require 
further stu dy of the paramete r space before they can be formulated. 

From lBate et al.l ( l2003h . the migration regime depends on the 
embryo mass. Type I migration occurs if M ^ Aft, and Type II if 
M > Aft, where the transition mass 



Mt = 2Af. ( ^ 



For Type I migration, the migration timescale is 

ti{a) 



M 



M, 



n 



a 



(45) 



(46) 



and for Type II 



(47) 



(where we have assumed a self-gravitating, Q = 1 disc). As ex- 
pected, type I migration is typically several orders of magnitude 
faster than type II. 

Planet-disc interactions do more than modify the planet's 
semi-major axis. Other orbital elements can also be modified, prin- 
cipally the eccentricity and inclination. Co-orbital Lindblad res- 
onances will typical ly act to dampen eccentricities if the planet 
is fully embedded t Artvmowica 1 19931) . but if the planet opens a 
gap dampening is suppressed, and can even result in eccentricity 
growth d epending on the strength and saturation of corotation res- 
onances JOgilvie &Lub ow 2003; Moorhead & Adams 2008) . In 
the case of massive planets (i.e. above I A fjup), simulations in- 
dicate that eccentricit y is generally damped dBitsch & Klevll20ICI ; 
I Kiev & Nelson! 1 20 I2I ). but a multi-embryo environment may pro- 
duce larger eccentricities through clu mp-clump scattering, or even 
eject them from the system entirely (Bolev et al.l20I0l ; IStamateIlosl 
l20I3b . We assume that the embryos adopt coplanar, circular orbits 
for simplicity. 



2.1.9 Tidal disruption 

If the embryo's radius exceeds its Hill Radius 

/ Af ^ 1/^ 



Rh 



\3AL 



(48) 



the upper layers can become unbound. As the semi-major axis de- 
creases, the radial acceleration will typically increase, while the 
embryo contraction rate decrease s with time . The in evitable result 
is tidal disruption. Appendix A of lNavakshinI ( 1201 lah indicates that 
while irradiation of the upper envelope by the star can also achieve 
disruption, tidal disruption tends to occur first for embryos around 
low mass stars (although the modification of the disc structure due 
to the radiation fie ld, not modelled in this work, can affect the mi- 
gration process, cf lFouchet et al.l20I2l) . 

The fate of the solid material inside the embryo depends on the 
evolutionary stage. As previously mentioned, if the grain sphere is 
not self-bound, then it will become unbound as the envelope disap- 
pears. If the embryo is completely dissolved, the resulting solids, 
which may have grown significantly since the embryo's birth, are 
disperse d into the protoplanetary d isc, as processed solids and rub- 
ble piles iBolev et alj feOlOl . uOllh . potentially even planetesimals 
JNavakshin & Chall2012l) . 

If the grain sphere has become self-bound and initiated its final 
collapse, then the stripping of the envelope can allow the emergence 
of a terrestrial mass protoplanetary core. 



2.2 Evolution of the Protostellar Disc 

We assume that the disc is axisymmetric and that the surface den- 
sity, 'E{r,t), evolves viscously according to equation ([3}. In a 
steady-state, this can be integrated to give an expression for the 
mass accretion rat e, Af , which, at radii large compared to the ra- 
dius of the star, is ( |Pringlell98Ir) 



Af = 37r!/E. 



(49) 



One of the main uncertainties in disc evolution is what pro- 
vides the kinematic viscosity. As discussed earlier, we assume that 
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the v iscosity has the form of an a viscosity JShakura & Sunvaevi 
1 1973) so that u = acsH, where Cs is the disc sound speed, H is 
the disc scaleheight (H = Cs/Q, with il the angular frequency of 
the disc), and a << 1 is a parameter that determines the efficiency 
of the angular momentum transport. 

We ran two sets of 100 disc models, selecting the central star 
mass randomly between M« = 0.8 M0 and A/* = 1.2 Mq. We 
use a ratioed grid that extends from r = 0.1 au to r = 1000 au. 

In one set of models, we assume the surface density of the disc 
model is initially non-zero inside r = 50 au only. When used in the 
population synthesis model, we assume this disc initially extends 
to some value of r randomly selected between 50 and 100 au, and 
that the fragmentation process truncates the disc to 50 au. 

In the second set of models, the disc initially has non-zero 
surface density inside its full 100 au extent, i.e. we assume that the 
fragmentation process does not truncate the disc. 

We assume an initial surface density profile of S oc r^^. Ide- 
ally, there should be a range of profiles measured: making the pro- 
file shallower may favour fragmentation by placing more mass at 
grea ter distance, which ty pically allows more rapid radiative cool- 
ing JForgan & Ricell2012h . In this analysis, we decided to choose 
a single profile to keep the scope of the work under control, but 
future work should investigate this. 

In each simulation, the initial disc mass is randomly chosen to 
be between 0.125A/, and 0.375A/*. The initial disc mass is there- 
fore quite high and such discs are likely to be self-gravitating. Note 
that we do not add mass to the disc in any way - this is equivalent to 
the implicit assumption that the initial phase of cloud collapse and 
protostar formation is complete, and that the accretion of material 
onto the disc is much less than the mass loss due to winds. This 
is an important assumption, as it curtails the accretion rates of the 
embryos, and future modelling work should address this. 

Gravitationally unstable discs can either fragment to form 
bound objects JBossll 19981 1200 2) or they can settle into a quasi- 
steady st ate in which the instability acts to transport angular mo- 
mentum JLin & Pringlel 19871 : iLaughlin & Bodenheimej|l994 . 

As made clear previously, it is now reasonably well estab- 
lished that fragmentation is unlikely in protostellar discs within 
20 — 30 au. If fragmentation does not take place, the gravita- 
tional instability is likely, at early times at least, to transport an- 
gular momentum outwards, allowing mass to accrete onto the cen- 
tral star. In such a state, the disc will achieve thermal equilibrium 
with dissipation due to the gravitational instability balanced by 
radiative cooling. Assuming Q = 1.5 therefore fixes the sound 
speed, and hence temperature, profile in the disc. I t has been shown 
JBalbus & Papaloizoiilll999l : lLodato & Ricell2004) that the gravita- 
tional instability transports angular momentum in a manner anal- 
ogous to viscous transport. The rate at which viscosity generates 
dissipation, D(R), in the disc 



D(R) = -uT.n^ 
4 



is balanced by radiative cooling. A, where JHubenvll 19901) 



A 



16 



-a T 



°> 1 



+ t2 



(50) 



(51) 



In Equation ([STJ, a is the Stefan-Boltzmann constant, Tc is 
the midplane temperature determined from Q = 1.5, To is a tem- 
perature due to external irradiation, and r is the optical depth. We 
assume r = kE, wh ere k, is the opacit y, and we use the analytic 
approximations fr om Bell & Lirj ( Il994f). As described in detail i n 
iRice & Armitag3 ( |2009|) (see also lciarkell2009l: Izhu et al.ll2009l) . 



the above can be used to determine the effective value of a. Equa- 
tion l3j can therefore be integrated to determine the time evolution 
ofS. 

If, however, the gravitational ct is less than 0.005, we assume 
that another trans port mechanism, such a s the magnetorotational 
instability (MRI) JBalbus & Hawle\ilT99ll) will then dominate and 
we set a = 0.005. In this case, we no longer require that Q = 1.5 
but set the midplane temperature, Tc, such that the cooling rate 
given by Equation \5\\ matches the dissipation rate given by Equa- 
tion l lSOb . In addition, we also assume that irradiation from the cen- 
tral star sets a minimum temperature in the disc. 

We also assume that radiation from the central star results 
in photoevaporative mass-loss, producing a disc wind (Ewind). 
There are a number of models for photoevaporative mass loss, 
many of which assume tha t UV photons from the central star 
ionise the disc surface (e.g ., iHollenbach et alJI 19941 : Iciarke et al 
l200l|: l Alexander et al.l2006n . Recent models, however JOwen et al 



2010), suggest that x-rays are the dominant driver and so here we 
implement the x-ray photoionization model described in detail in 
Owen etal. (2011). 

Each disc model uses a different star mass, but assumes the 
same disc-to-star mass ratio. Despite this, the variation in x-ray 
luminosity (from 5 x lO'^* erg s~^ to 10^^ erg s~^) produces a 
wide range of different disc lifetimes, consistent with that observed 
(Haisch etal. 2001). 

We therefore have a disc model that can self-consistently 
evolve the surface density from the early stages when the gravi- 
tational instability is likely to dominate through to the later stages 
when an alternative transport mechanism, such as MRI, will dom- 
inate and also includes the late-stage dispersal due to photoevapo- 
rative mass-loss. 



2.3 Implementation 

To conclude this section, we summarise the algorithm of the semi- 
analytic model. 

For each planetary system, we begin by randomly selecting a 
disc model. We either assume the disc truncates after fragmenta- 
tion, or does not. If we are dealing with a truncated disc model, we 
initially extend it to a maximum radius rmax beyond the model's 
outer radius of 50 au. We sample rmax uniformly between 50 and 
100 au, and we assume the disc surface density profile remains con- 
stant in this extension, and calculate the subsequent disc proper- 
ties accordingly. If we are dealing with a non-truncated disc model 
(which extends to 100 au), we do not extend it in any way. Again, 
rmax is somewhat arbitrary - to first order, increasing rmax will in- 
crease the number of fragments initially. Whether this will increase 
the number of fragments that survive the tidal downsizing process 
is unclear. Future work should investigate how changing this affects 
the resulting distribution of objects. 

The Jeans criterion is then used to find the minimum disc ra- 
dius at which fragmentation is expected to occur. Embryos are then 
created, with mass equal to the local Jeans mass, at spacings of a 
few Hill radii (where the exact value is randomly sampled). The 
other initial embryo properties are calculated according to the first 
core equations. If we are assuming disc truncation, then we revert 
back to the initial model, removing the material added to extend it 
to rmax. If we do not truncate the disc, then we make no modifica- 
tions at this stage. 

The entire planetary system is then evolved on a timestep con- 
strained by the physical timescales in play. The timestep can be no 
larger than the minimum value of fgr, fsod or icooi for any embiyos 
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in the simulation. The disc model is interpolated from a look-up ta- 
ble to update its properties with time, and these are used to calculate 
the migration timescales of all embryos. Note that each embryo is 
evolved independently of each other - while several embryos may 
be evolved in the same disc in parallel, they do not interact in any 
way, an important limitation (see Discussion). 

Once an embryo's igr (which will be a function of Scrit, itself 
a function of disc radius) has elapsed, sedimentation is allowed to 
proceed, and the size of the grain sphere Rs is evolved (with the 
sedimentation velocity not allowed to exceed 10 metres per second 
due to fragmentation of the grains). 

At every timestep, the grain sphere is tested to check if 
it is self-gravitating, or self-bound. Grain spheres that are self- 
gravitating use the accelerated sedimentation equation, and grain 
spheres smaller than the Jeans length are considered to be self- 
bound, and immediately form solid cores. 

Throughout this grain growth process, the embryo's central 
temperature is monitored. If the embryo exceeds temperatures of 
180 K before core formation, the ice mantles are melted, and two 
thirds of the solids mass is deducted from the embryo. If the tem- 
perature exceeds 1400 K, the dust is vapourised and core formation 
is completely halted. Embryos with temperature above 2000 K en- 
ter the second core phase and are assigned final radii according to 
the mass-radius relation of burrows et al.l ( 119971) . 

As the embryo migrates inward, it will begin the process of 
tidal disruption. Instead of immediately stripping the embiyo, the 
embryo's upper layers are removed more gradually over the course 
of one orbital period, beginning at time tstrip with the radius de- 
caying exponentially to its Hill radius Rh according to: 



R{t) = R{ts 



exp 



{Rit: 



strip/ 



Rh) 



(t — tstrip) " 
2^ 



(52) 



The mass M of the embryo is then updated assuming the embryo 
remains a polytrope of index 3/2. If the new embryo radius is less 
than the grain sphere radius, then the total mass in solids is also 
reduced by an appropriate amount. 



3 RESULTS 

To investigate the dependency on initial input parameters, we run 
the population synthesis model with several different configura- 
tions: 

(i) Setting Pk = 1, and using disc models which initially extend 
to a radius randomly selected between 50 and 100 au, and are then 
immediately truncated to 50 au after fragmentation has occurred. 

(ii) As 1., except that the migration timescales are uniformly in- 
creased by a factor of 10, i.e. migration is ten times slower 

(iii) As 1., except p^ = 2. 

(iv) As 1., except the disc models extend to 100 au, and are not 
truncated after fragmentation. 

Every disc system is run for at least IMyr, to ensure the full span 
of disc and embryo evolution is observed. 



3.1 Opacity Law Index p^ ~ 1 (Truncated Disc) 

The left panel of Figure [T] shows the mass-semimajor axis distri- 
bution of the embryos at the point of fragmentation. As expected, 
the embryos initially occupy mass ranges above a few Afjup, and 



fragmentation does not occur at small a. As some of the discs are 
quite massive and cool, the fragmentation radius can be as low as 
20 au. The colour bar indicates the grain sedimentation timescale 
for the embryos, which in the majority of cases is relatively short 
(< 200 yr). If the grains can grow sufficiently quickly that sedi- 
mentation can begin before vapourisation and disruption, core for- 
mation seems possible for many of the embryos. 

However, the rapid temperature evolution of the embryos, cou- 
pled with tidal disruption beginning almost immediately after frag- 
mentation, tends to prevent core formation. The right panel of Fig- 
ure[T]shows the mass-semimajor axis distribution of the fragments 
foiTned in the population synthesis model after 1 Myr. Note that 
most embryos do not possess a solid core, and are sufficiently mas- 
sive to be considered brown dwarfs. The disc inner boundary is 
located at 0.1 au, and as such all objects which reach this location 
without being tidally disrupted are not evolved any further. 

There is a lack of objects at semimajor axes below ~ 
5 au - somewhat consistent with the "brow n dwarf desert" 
JMarcv & Bufleil I2OO0I: iHalbwachs eTal] l2000t) ) - which can be 
seen in the right hand panel of Figure |2] Generally, objects above 
the brown dwarf mass limit exist at semimajor axes > 20 au, 
showing that in situ formation of brown dwarfs via fragmenta- 
tion remains viable in the tidal do wnsizing paradigm, in agreement 
with recent numerical simulations (IStamatellos & WhitworthI 20091 : 
Kaplan etal. 2012 ). 

If we consider both histograms in Figure |2] we can see that at- 
trition of the embryo thanks to tidal disruption, and the loss of em- 
bryos through migration strongly shapes the final statistics of the 
objects. Around 40% of all embryos formed are lost - as the mass 
of the embiyo decreases during the disruption process, the migra- 
tion regime switches from Type II to Type I, significantly reducing 
the migration timescale and forcing objects towards the central star. 
The mass distribution of objects (left panel of Figure |2j shows a 
significant number of objects formed above SOAfjup, with masses 
in the hydrogen-burning M dwarf regime, but the final mass distri- 
bution tends to fall off above this value as disruption acts on these 
embryos before nuclear reactions begin. 

While there are some embryos with final semimajor axes be- 
low 10 au (right panel of Figure|2), the relative frequency of objects 
is quite low. Indeed, the final semimajor axis distribution retains 
a strong signature of the minimum fragmentation radius of about 
20 au. Despite this fragmentation boundary being much closer to 
the central star than the values of 40 — 70 au typically found in 
numerical simulations, the population of close-in objects remains 
quite small compared to the large semimajor axis population. This 
appears to be due to disc truncation - as most of the objects will 
form in regions that have experienced truncation, they will not be- 
gin migration until the disc has spread on viscous timescales. This 
increases the migration timescale, and allows them to evolve with 
limited tidal disruption and reach the "second collapse" due to H2 
dissociation before the disc can act to move the objects inward, and 
increase the rate of tidal disruption. 

The left panel of Figure [3] shows the distribution of solid core 
masses formed in this run. The distribution is peaked at around 
6M(g, with a tail extending as far as 22M^. The extrema of the 
distribution compares well to (fo r example) the range of estimate s 
for Jupiter's core mass (0-18 Mm . iGuillotl 1 9991 : lNettelmannl20 1 lb . 
but we must also remember that objects formed in the model with 
these core masses are likely to have total masses of lOMjup or 
above. 

The right panel of Figure |3] illustrates that brown dwarfs are 
the most populous of all the surviving objects in this run (around 
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Figure 1. Embryo mass versus semi-major axis for p^ = 1, where tlie disc truncates at 50 AU after fragmentation. The left hand panel shows the initial 
fragment mass-semimajor axis distribution, with the colour bar indicating the timescale on which grains will sediment once sufficiently grown. 256,424 disc 
fragments were produced in this mn. The right hand panel shows the final distribution after lO'' years. Of the fragments produced, only 116,150 survived the 
tidal downsizing process, and 50,543 were able to fonn solid cores. The colours in this plot indicate the mass of the core in Earth masses. The disc model's 
inner boundary is fixed at 0. 1 au, and as such objects cannot cross this boundary. 
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Figure 2. Comparing the initial (blue lines) and final (red bars) distributions of embryo mass (left) and embryo semimajor axis (right). 



50% of the total), closely followed by objects with a solid core (at 
around 42%). No icy cores are formed, as the embryos quickly ex- 
ceed the melting temperature before core formation. If and when 
the grain growth process is successful, the sedimentation timescale 
is sufficiently short that cores are formed before full disruption oc- 
curs. We should therefore conclude that (for these parameters) these 
disruption events do not add significant processed materials or sub 
Earth mass solid bodies to the surrounding environment. 



3.2 Opacity Law Index Pk = 1 (Truncated Disc, Inefficient 
Migration) 

In the same vein as llda & Linl (2003), ws now vary the speed of 
migration by multiplying the migration timescales in equations l l46t 
and l |47t by a factor of 10. Figure |4]shows the initial and final mass 
versus semimajor axis distributions in this run. 

With migration ten times less efficient than the previous run, 
semi-major axes lower than the fragmentation boundary are less 
populated than before. Some objects still migrate to the disc inner 
boundary, but in much fewer numbers than before. The sedimen- 
tation timescale for the fragments remains the same, so the cores 



produced are very similar to the previous case. With the ratio be- 
tween sedimentation and migration timescale being reduced by a 
factor of ten, the formation of solid cores before disruption is some- 
what easier. Indeed, this mn shows the only example of a terrestrial 
planet being formed in the simulation, with mass less than 1 M^ at 
a semimajor axis of 1 au. 

The embryo mass distribution (Figure |5] left panel) remains 
relatively unchanged, but the semimajor axis distribution (right 
panel) shows an even larger desert, with few objects ending up in- 
side 20 au. 

The core mass distribution (left panel of Figure|6) is again very 
similar as the sedimentation behaviour is unchanged. Brown dwarfs 
still dominate the population (right panel of Figure |6), followed 
closely by cored objects. 



3.3 Opacity Law Index p^, — 2 (Truncated Disc) 

In this run, we now allow migration to be efficient, and truncate 
the disc at fragmentation, but change the opacity law so that it 
scales with the square of the temperature. This significantly in- 
creases the cooling time of the embryo, which will have important 
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Figure 3. Left: The distribution of core masses after 10 years of evolution. Right: The relative frequency of each object type formed in the population 
synthesis model, (light). 





Figure 4. Embryo mass versus semi-major axis for p^ = 1, vi'ith a disc that truncates at 50 AU after fragmentation and migration timescales increased by a 
factor of 10. The left hand panel shows the initial fragment mass-semimajor axis distribution, with the colour bar indicating the timescale on which grains will 
sediment once sufficiently grown. 202,385 disc fragments were produced in this run. The right hand panel shows the final distribution after 10® years. Of the 
fragments produced, only 92,295 survived the tidal downsizing process, and 40,571 were able to form solid cores. The colours in this plot indicate the mass of 
the core in Earth masses. 
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Figure 5. Comparing the initial (blue lines) and final (red bars) distributions of embryo mass (left) and embryo semimajor axis (right). 
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Figure 6. Left: The distribution of core masses after 10 years of evolution. Right: The relative frequency of each object type formed in the population 
synthesis model, (light). 
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Figure 7. Embryo mass versus semi-major axis for p^ = 2, with a disc that tmncates at 50 AU after fragmentation. The left hand panel shows the initial 
fragment mass-semimajor axis distribution, with the colour bar indicating the timescale on which grains will sediment once sufficiently grown. 530,488 disc 
fragments were produced in this run. The right hand panel shows the final distribution after 10® years. Of the fragments produced, only 22,637 survived 
the tidal downsizing process, and 2,292 were able to form solid cores. The colours in this plot indicate the mass of the core in Earth masses. The left hand 
panel shows the initial fragment mass-semimajor axis distribution, with the colour bar indicating the timescale on which grains will sediment once sufficiently 
grown. 



consequences for the population produced. Figure |7] shows the ini- 
tial and final mass vs semimajor axis distributions for this run. 

The inefficient cooling of the embryo results in a greatly in- 
creased destruction rate - over 95% of the fragments produced in 
this run are destroyed. Only the most massive objects survive, and 
there are only a handful that survive in the inner disc. This is re- 
flected in the mass and semimajor axis distributions (Figure[8]l. 

Perhaps most worrying are the extremely large solid core 
masses produced in this run (left panel of Figure |9). The mode of 
this distribution is around 200i\f® ~ O.SAfjup. Note that these 
very large cores only exist in massive gaseous embryos - there 
are no terrestrial type objects (i.e. solid cores without gaseous en- 
velopes). These anomalously large cores are suggestive of missing 
physics from the model (see Discussion). 

The inefficiency of embryo cooling only allows massive ob- 
jects to remain bound without becoming disrupted, and as a result 
cores form very rarely in this run (only around 10% of all surviv- 
ing objects). The other 90% form brown dwarfs - indeed, some 
are massive enough to form hydrogen-burning low-mass stars with 
masses ~ O.IMq. 



3.4 Opacity Law Index Pk ~ i (No Disc Truncation) 

We have seen in the previous runs that the migration of the embryos 
strongly determines their mass and semimajor axis evolution. Until 
now, we have assumed that the fragmentation process truncates the 
disc, which allows embryos to evolve their physical radius and dust 
population before the disc can undergo viscous spreading, extend- 
ing its outer radius until it encroaches upon their orbital domains, 
and begin the migration process. 

In this run, we do not truncate the disc at fragmentation. The 
disc can therefore begin migrating the embryos immediately, de- 
creasing the Hill Radius of the embryos precipitately and signifi- 
cantly increasing the rate of disruption. As the p^ = 2 run already 
shows a high destruction rate, we reset p^ to 1 to assess the signif- 
icance of disc truncation. 

Figure [To] shows the initial and final mass vs semimajor axis 
distributions, which indicate that removing disc truncation has a 
significant effect. Again, there are no terrestrial type objects, and 
all objects which possess a core possess a massive envelope. With 
migration occuring immediately after fragmentation, many more 
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Figure 8. Comparing the initial (blue lines) and final (red bars) distributions of embryo mass (left) and embryo semimajor axis (right). 
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Figure 9. Left: The distribution of core masses after 10*' years of evolution. Right: The relative frequency of each object type formed in the population 
synthesis model, (light). 
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Figure 10. Embryo mass versus semi-major axis for p^ = 1, with a disc that does not truncate after fragmentation. The left hand panel shows the initial 
fragment mass-semimajor axis distribution, with the colour bar indicating the timescale on which grains will sediment once sufficiently grown. 205,184 disc 
fragments were produced in this run. The right hand panel shows the final distribution after 10® years. Of the fragments produced, only 74,430 survived 
the tidal downsizing process, and 7,378 were able to form solid cores. The colours in this plot indicate the mass of the core in Earth masses. The left hand 
panel shows the initial fragment mass-semimajor axis distribution, with the colour bar indicating the timescale on which grains will sediment once sufficiently 
grown. 
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objects reach the inner disc boundary, producing a small population 
of Hot Jupiters and close BD companions. Core formation is sup- 
pressed compared to the other p^ = 1 runs: cores form in 10% of 
surviving embryos compared to around 40% in the other p^ cases. 

The mass distribution of the objects (left panel of Figure [TT) 
remains relatively unchanged from initial to final, forming large 
numbers of objects up to 100 A/jup. The two modes of the distri- 
butions are generated by the spacing of the fragments as they form 
in the disc (with some broadening of the peaks due to the uniform 
sampling of the spacing in Hill Radii described in earlier sections). 

This population is very much a stellar population. The brown 
dwarf desert is less significant (right panel of Figure 1111 . and few 
objects survive with semimajor axis larger than 60 au. 

The core mass distribution (left panel of Figure 112) also re- 
flects the initial mass distribution. Less than 10% of all surviving 
objects have cores (right panel of Figure [Tzt. and nearly 90% of all 
the objects are brown dwarfs (or low mass stars). 



4 DISCUSSION 

4.1 Limitations of the Model 

We readily admit that this model, while detailed in construction, 
still lacks some important physics, which we hope to include in 
future work. 



dNavakshinluOl id) . This would suggest that embryos rotate suffi- 
ciently rapidly to form c ircum-embryonic discs (Bole y et al. 2010|; 
IShabram & BolevI 120131) and binary cores, an outcome not pro- 
ducible in this current model. 

Equally, the rotation rate of giant planets is expected, as a rule, 
to be initially prograde, and close to the break-up frequency, al- 
though stellar tides can damp this rotation rate. Also, events such 
as planet-planet scatt ering or the Kozai me chanism can produce 
retrograde orbits (e.g. iNagasawa & Idall201 ih . Future work should 
model this angular momentvrai evolution in detail. 



4.1.2 Dust Coagulation 

The extremely massive cores produced by the p^ = "2 run sug- 
gest that there is missing physics in the coagulation of the dust 
grains. While turbulence is part of the coagulation process in these 
models, convection is not. As we have mentioned previously, al- 
lowing the envelope to b e convective can sig nificantly increase the 
sedimentation timescale JHelled et alJuOOSi) , although it does not 
completely preclude core formation, especially for low-mass em- 
bryos (Helled & Schubert 2008). Also, we should strictly alter the 
opacity in the embryo as the grains sediment, as the local metal- 
licity of the gas can change significantly. Global metallicity effects 
are also important: decreasing the initial dust to gas ratio in the 
disc can greatly decrease the contraction timescale of the embryo 
(iHelled & Bodenheimern2011i) . 



4.1.1 Accretion of Mass and Angular Momentum by the Embryo 

Perhaps the most glaring omission from this population synthesis 
model is accretion of material onto the embryos after fragmenta- 
tion. As the semi analytic disc model is evolved without the pres- 
ence of fragments, modelling subsequent mass loss from the disc 
onto embryos is not possible. The embryo should be able to ac- 
crete provided it does not open a significant gap, which it will do 
if a) its Hill radius exceeds the local disc scale height, and b) the 
clump is sufficiently massive to overcome the local Reynolds turbu- 
lence JLin & Papaloizoul 19791) . For almost all fragments formed, 
this criterion is easily satisfied initially. However, as the embryo 
migrates inward, the Hill radius decreases, and the embryo's mass 
also decreases due to tidal disruption. The gap opening criterion 
will no longer be satisfied, and accretion can begin once more. 

To what level this latter accretion phase will affect the subse- 
quent evo l ution of fragments is unclear. Numerical simulations by 
Izhu et al.l ( 12012 ) show gap opening at semimajor axes as low as 10 
au, which would be a consequence of accretion increasing the Hill 
radius past the value of the local scale height. 

We have also ig nored the embryo' s ability to accrete solids 
from its surroundings (IHelled et al .120061) . The embryo may be able 
to substantially enrich itself during the early stages of its contrac- 
tion, accreting a large fraction of the available planetesimals in its 
feeding zone, boosting the dust-to-gas ratio. This could be incorpo- 
rated into our model quite simply using simple two-body capture 
approximations. 

The ac cretion history of the emb ryo will determine its angular 
momentum. ICha & NavakshinI ( 120111) find in their numerical sim- 
ulations that the fragments produced rotate prograde to the disc's 
rotation, with a rotational angular frequency approximately 10% of 
their maximum break-up frequency. This can be used to estimate 
the rotational angular frequency of terrestrial planets formed from 
these embryos whose maximum exceeds the break-up frequency 



4.1.3 Embryo Migration 

We have deliberately implemented a basic model of planetary mi- 
gration in this work, as this allows us to investigate the com- 
plete system of equations presented by Navakshin (2010b, 201 la, bj) 
and assess their expected output on a statistical level. As such, 
the modelling of migration misses out important characteristics 
of the migration pro cess: stochasticity due to disc turbulence 
(Baruteau et al. 2011), resonant migration of multiple bodies (e.g. 
iLibert & Tsigani s 201 1) , an d halting of migratio n at planet traps in- 
duced by disc structure (e.g lKretke & Linll2012h . Inwardly migrat- 
ing giant planets can even halt or reverse the motion of low mass 
planets at larger semimajor axis dPodlewska-Gaca et al.ll2012h . 

Also, the eccentricity and inclination evolution of the bodies is 
not tracked: recent numerical simulations suggest that eccentricity 
excitation via planet-dis c interaction will be low for planet masses 
of a few Jupiter masses dDunhill et al. 20121). w hich appears to be 
borne out by observations (Dawson et al.ll2012h . but future models 
should investigate this by allowing eccentricities to be excited and 
damped. 

Most migration models are developed for low-mass, non- 
self- gravitating discs: simulations of fragment migration in self- 
gravitating discs have shown that the usual embryo migration 
equations do not appe ar to apply, and that fragment migration is 
typically much faster ( iMichael et"al]|201ll : iBaruteau etai]|201ll : 
IZhu et alj|2012h . Also, as our embryos are evolved independently 
of each other, we cannot model for example the chaotic orbital in- 
terchanges wit nessed in se lf-gravitating disc simulations with mul- 
tiple embryos ( lBossll2013h . It is unclear how these interactions can 
be modelled semi-analytically, but it is clear that any model output 
should be considered with the absence of this physical process at 
the forefront of one's mind. 

Can t hese simulations produce semi-analytic prescriptions for 
migration? IBaruteau et al.l ( 120111) show that gravitoturbulence re- 
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Figure 11. Comparing the initial (blue lines) and final (red bars) distributions of embryo mass (left) and embryo semimajor axis (right). 
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Figure 12. Left: The distribution of core masses after 10 years of evolution. Right: The relative frequency of each object type formed in the population 
synthesis model, (light). 



duces the Type I migration timescale by a factor of a few as a re- 
sult of the disc's non-axisymm etry pushing the lo cation of Lindblad 
resonances towards the planet. IZhu et al.l (120 12h also show that the 
Type I timescale is altered by a factor which they fit in terms of em- 
bryo mass and local surface density, accounting for the embryo's 
inertia. Future modelling should consider implementing some of 
these initial empirical results. 



4.2 Terrestrial Planet Formation via Tidal Downsizing 

Despite the above limitations, we can make some reasonably robust 
statements regarding the likelihood of terrestrial planet formation. 
All the runs carried out in this paper show that planets of mass 
~ Afg are extremely rare - out of over a million disc fragments, 
only one stands out as a clear example of a terrestrial planet. 

Most of the missing physics described above is likely to inhibit 
the formation of low mass cores rather than encourage it, so future 
population synthesis modelling is likely to produce even fewer ter- 
restrial planets than we see here. In short, it appears to be clear that 
while GI followed by tidal downsizing can produce terrestrial plan- 
ets on some occasions, it is not an efficient means of constructing 
such objects. This would indicate that terrestrial planets are much 
more likely to be formed by core accretion than disc fragmentation. 



4.3 Production of Asteroid Belts 

For the production of asteroid belts via tidal downsizing, the grain 
growth process must be well advanced, without the sedimentation 
process forming a solid core or the embryo temperature increasing 
above the vapourisation temperature. If the disruption process is 
sufficiently rapid, then the embryo and core will unbind, and the 
grains in the embryo will be distributed in a ring at the semimajor 
axis where disruption began (Nayakshin & Cha 2012). 

In these models, we did not see any embryos which com- 
pletely satisfy the above conditions for belt formation. Only one run 
(the non-truncated p^ = 1 run) showed some embryos with grains 
surviving the disruption process without either being vapourised or 
forming a core. The embryos showed evidence of ice melting, but 
their grains had not grown to the point where sedimentation would 
begin. As such, these would not form asteroid belts in the conven- 
tional sense, but could perhaps for m dust belts with some grains 
already grown (cf lBolev et al.ll201(]|) . 

This would skew the disc grain size distribution towards larger 
values, potentially providing a "shortcut" to grain growth in the 
core accretion model. This enhanced grain growth is similar to the 
grain growth observed in non-fragmenting self-gravitating disc spi- 
ral structures, which has also been proposed as a means of accel- 
erating the core accretion process (IRice et al.ll2004t iGibbons et al.l 
i2012) . The key difference between embryo-enhanced grain pop- 
ulations and spiral-wave enhanced populations is the location of 
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the grains' deposition: spiral-wave enhianced grains will occupy 
semi major axes above several tens of au ( Clarke & Lodato 20093, 
whereas embryo-enhanced grains can be delivered to lower semi 
major axis before the embryo is destroyed. 

If sedimentation has just begun when the embryo is totally 
disrupted, then a large fraction of the grains supplied to the disc 
may be of order millimetre si zes, allowing rapid core accretion 
dLambrechts & Johansenll2012h . 



4.4 Embryo disruption and Protostellar Outburst Models 

In this work, typically around 40% of all embryos were com- 
pletely destroyed by tidal disruption. Thi s large destruction rate , 
borne out by num er ical simulations (e.g. J Bolev & DurisenI 120101 : 
lNavakshirJ|2010Jlbl [JoTTal : Izhu et alj|2012h . could result in pro- 
tostellar outburst behaviour such as the FU Orionis phen omenon 
jVorobvov & Basul2005ir2006l : lDunham & Vorobvovll2012 ). If this 
is the case, then the behaviour of the gas as it becomes unbound 
from the embryo will determine the nature of the burst. Indeed, 
subsequent disc fragmentation can be suppressed by the star's out- 
burst (Stamatellos et al. 2012). We have not coupled disc evolution 
and embryo evolution, and so we cannot model such behaviour. As 
a result, the population synthesis model does not characterise the 
outbursts, but merely suggests that episodic stellar accretion is a 
common outcome of the disc fragmentation process. 

More subtly, Navakshin & Lodato (2012) show how a tidally 
disrupting embryo can create deep gaps in the inner disc that are re- 
filled by its own matter, which can then be accreted onto the central 
star, producing outbursts and quiescent p eriods as the gap refills , 
without the embryo's destruction (see also lLodato & Clarke 2004). 
While this population synthesis model does not self-consistently 
evolve the embryo evolution and the disc evolution, this phe- 
nomenon could be modelled semi-analytically in the future. 



4.5 Brown Dwarfs and the Desert 

The population synthesis model predominantly produces brown 
dwarfs - this is true for all four runs carried out in this work. De- 
pending on the opacity law, as much as 90% of the objects that 
survive the tidal downsizing process satisfy the canonical M > 
ISAfjup criterion for being brown dwarfs. Doe s the brown dwarf 
desert (iMarcv & Butler 2000: Halbwachs et al.i 2000) survive the 
migration process? In Figure [13] we plot the semimajor axis distri- 
bution of the brown dwarfs formed in all four runs. 

If we define the desert as a lack of objects below around 5 au, 
then the runs with p^ = 1 and disc truncation show this lack quite 
strikingly. If migration is suppressed (top right of Figure [T3ll, then 
this desert extends as far as the fragmentation boundary of the disc 
model. Changing the opacity law or not truncating the disc makes 
the desert less striking - brown dwarfs are then able to survive to 
the inner edge of the disc model at 0.1 au. The Pk ~ 2 run still 
preserves a desert in between 0.2 < r < 20 au, but the run without 
disc truncation shows a desert with greatly reduced significance. 

This would suggest that for a strong brown dwarf desert to 
be formed, the disc should truncate during the fragmentation pro- 
cess (or migration should be inefficient), and the embryo's opacity 
should not increase too quickly with temperature. We do not fully 
model the so-called "opacity gap" normally se en at temperatures 
between ~ 500 -5000 K (e.g. Bell & Lirj[T994l ). which would sig- 
nificantly affect the collapse of the embryo envelope. A new set of 
grain growth and sedimentation equations where p^ is allowed to 



vary with temperature would also significantly affect the potential 
for core formation. 

Also, we do not model ejection o f brown dwarfs 
via scattering or three-bod y interactions (iBate et al.l l2003l : 
I Goodwin & WhitworthI 120071) . which may be more e fficient at 
removing coreless objects than objects with cores (ILiu et al.l 
I2OI2I) . Regardless, our results would suggest that the destruction 
of inwardly migrating brown-dwarf-mass objects is already quite 
efficient at remo ving close separation brown dwarf companions 
JArmitage & Bo nnell 2002), and that the restriction of disc frag- 
mentation to semimajor axes well above 5 au would ensure those 
brown dwarfs which survive disruption would be unlikely to form 
close separation binaries with their parent star. 



5 CONCLUSIONS 

We have presented the first attempt at a population synthesis model 
for objects produced by self-gravitating disc fragmentation and 
subsequent tidal downsizing. The model imposes a protoplanetary 
disc, which can then fragment under gravitational instability. The 
fragments can then migrate, grow grains, sediment those grains to 
the centre of the fragment and form solid cores, or be tidally dis- 
rupted (not necessarily in that order). Each protoplanetary disc (if 
it fragments) produces several objects - running the model many 
times with different disc properties produces a large number of ob- 
jects, allowing the construction of statistics on tidal downsizing as 
a planet formation mode. 

We ran the population synthesis model for four different sce- 
narios, varying the speed of migration, how the opacity within the 
embryo changes with temperature, and whether the fragmentation 
process truncates the disc or not. Each run simulates a sufficient 
number of discs to produce ~ 10^ fragments. 

In general, we find that around half of all fragments are com- 
pletely destroyed during the tidal downsizing process. The majority 
of surviving fragments remain at semimajor axes > 20 au after 1 
Myr of evolution. Their masses are typically greater than 5 i\fjup, 
with a large fraction exceeding the brown dwarf mass limit of 13 
Mjup- Solid cores form with a distribution of masses peaked typi- 
cally at 6 M(Q . Only a very small fraction of these solid cores lose 
their outer envelopes to form terrestrial planets (< 0.001%). We 
see no production of planetesimal belts via evolved grain popu- 
lations becoming unbound during the destruction of the embryo 
(although we saw some weak evidence of processed solids being 
returned to the disc via disruption). 

The model process does require further development, 
and the physics to be added will both suppress and en- 
hance grain growth/core formation. More appropriate mi- 
ration models may speed up the tidal disruption process 
.Baruteau et al.ll20Il[), but equally more fa ithful modelling of en- 
richment at birth ( IBolev & DurisenI I2OIOI) planetesimal capture 
(Helled et al. 200&) and opacity evolution during embryo contrac- 
tion (H elled & Bodenheimer 2011) can produce cores m uch more 
rapidly ( iHelled et al.l2008l : lHelled & Bodenheime3l201 ih . In short, 
it remains possible that the gravitational instability mode of planet 
formation can form a modest number of terrestrial planets, but it 
appears to form a much larger quantity of giant planets and brown 
dwarfs that preferentially exist at semimajor axes above around 30 
AU. 
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Figure 13. The final semimajor axis distribution of all brown dwarfs for the four runs carried out in this work. 
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